import numpy as np


class AlgorithmUtil:
    # 坐标转换（4326->3857）
    @staticmethod
    def WGS84ToMercator(WGS84_point):
        mercator_point = [None, None]
        earthRad = 6378137.0
        mercator_point[0] = WGS84_point[0] * np.pi / 180 * earthRad
        a = WGS84_point[1] * np.pi / 180
        mercator_point[1] = earthRad / 2 * np.log((1.0 + np.sin(a)) / (1.0 - np.sin(a)))
        return mercator_point

    # 坐标转换（3857->4326)
    @staticmethod
    def mercatorToWGS84(self, mercator_point):
        WGS84_point = [None, None]
        WGS84_point[0] = round(mercator_point[0] / 20037508.34 * 180, 10)
        y = mercator_point[1] / 20037508.34 * 180
        WGS84_point[1] = round(
            180 / np.pi * (2 * np.arctan(np.exp(y * np.pi / 180)) - np.pi / 2), 10
        )
        return WGS84_point

    # 计算余弦距离(向量夹角)
    @staticmethod
    def cosinDistance(vectorA, vectorB):
        # 向量点乘
        AB = vectorA[0] * vectorB[0] + vectorA[1] * vectorB[1]
        # 向量模长
        vectorA_moduleLength = np.sqrt(
            vectorA[0] * vectorA[0] + vectorA[1] * vectorA[1]
        )
        vectorB_moduleLength = np.sqrt(
            vectorB[0] * vectorB[0] + vectorB[1] * vectorB[1]
        )
        # 计算夹角余弦值
        cosin = AB / (vectorA_moduleLength * vectorB_moduleLength)
        if cosin > 1:
            cosin = 1
        if cosin < -1:
            cosin = -1
        angle = np.degrees(np.arccos(cosin))
        return 1 - cosin, angle
    # 计算轨迹向量序列
    @staticmethod
    def caculateVector(trajectory):
        # 投影计算
        projectedTrajectory = list(
            map(lambda point: AlgorithmUtil.WGS84ToMercator(point), trajectory)
        )
        # 初始化向量序列
        vectorList = []
        # 初始化向量角度序列
        angleList = []
        # 指北单位向量
        northPointingVector = np.array([0, 1])
        # 遍历每一个轨迹点
        for index in range(1, len(trajectory)):
            # 向量计算
            vector = [
                projectedTrajectory[index][0] - projectedTrajectory[index - 1][0],
                projectedTrajectory[index][1] - projectedTrajectory[index - 1][1],
            ]
            # 加入集合中
            vectorList.append(vector)
            vector = np.array(vector)
            # 计算向量夹角余弦值
            cos = vector.dot(northPointingVector) / np.sqrt(vector.dot(vector))
            # 判断当前方向向量所在象限
            if vector[0] < 0:
                # 向量在第二或第三象限
                angle = 360 - round(np.degrees(np.arccos(cos)), 3)
            else:
                # 向量在第一或第四象限
                angle = round(np.degrees(np.arccos(cos)), 3)
            angleList.append(angle)
        return vectorList, angleList
    # 计算轨迹的长度，并将轨迹坐标从4326转为3857，同时记录每两个轨迹点间的距离
    @staticmethod
    def getTrajectoryLen(trajectory):
        # 坐标投影
        projectedTrajectory = list(
            map(lambda point: AlgorithmUtil.WGS84ToMercator(point), trajectory)
        )
        # projectedTrajectory = trajectory
        # 初始化总长度和轨迹段长度列表
        totalLen = 0
        lenList = []
        # 遍历每一对坐标
        for index, point in enumerate(projectedTrajectory):
            if index == 0:
                continue
            # 计算两轨迹点间的长度
            curLen = np.linalg.norm(
                np.array(projectedTrajectory[index - 1])
                - np.array(projectedTrajectory[index])
            )
            lenList.append(curLen)
            # 长度累加
            totalLen = totalLen + curLen
        return projectedTrajectory, totalLen, lenList
